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This study presents an approximate model for the atypical Schumann resonance in Titan's atmosphere 
that accounts for the observations of electromagnetic waves and the measurements of atmospheric con- 
ductivity performed with the Huygens Atmospheric Structure and Permittivity, Wave and Altimetry 
(HASI-PWA) instrumentation during the descent of the Huygens Probe through Titan's atmosphere in 
January 2005. After many years of thorough analyses of the collected data, several arguments enable 
us to claim that the Extremely Low Frequency (ELF) wave observed at around 36 Hz displays all the char- 
acteristics of the second harmonic of a Schumann resonance. On Earth, this phenomenon is well known to 
be triggered by lightning activity. Given the lack of evidence of any thunderstorm activity on Titan, we 
proposed in early works a model based on an alternative powering mechanism involving the electric cur- 
rent sheets induced in Titan's ionosphere by the Saturn's magnetospheric plasma flow. The present study 
is a further step in improving the initial model and corroborating our preliminary assessments. We first 
develop an analytic theory of the guided modes that appear to be the most suitable for sustaining Schu- 
mann resonances in Titan's atmosphere. We then introduce the characteristics of the Huygens electric 
field measurements in the equations, in order to constrain the physical parameters of the resonating cav- 
ity. The latter is assumed to be made of different structures distributed between an upper boundary, pre- 
sumably made of a succession of thin ionized layers of stratospheric aerosols spread up to 150 km and a 
lower quasi-perfect conductive surface hidden beneath the non-conductive ground. The inner reflecting 
boundary is proposed to be a buried water-ammonia ocean lying at a likely depth of 55-80 km below a 
dielectric icy crust. Such estimate is found to comply with models suggesting that the internal heat could 
be transferred upwards by thermal conduction of the crust, while convective processes cannot be ruled 
out. 

© 2012 Elsevier Inc. All rights reserved. 


1. Introduction 

Titan, Saturn’s largest satellite, is the only moon in the Solar 
System with a dense atmosphere and stable liquid hydrocarbon 
lakes on its surface. Its global properties (mass, radius, density) 
are similar to those of its cousin icy jovian satellites Ganymede 
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and Callisto. The discovery by the Galileo mission of their own 
magnetic fields, partly intrinsically produced (Ganymede) and 
partially induced by the variations of the jovian magnetic field, 
has led to the interpretation that internal water oceans are present 
beneath their icy crust (Schubert et al., 2004). Discovering a buried 
ocean at Titan would mean that this moon is a likely world of 
organic molecules and liquid water, two of the necessary ingredi- 
ents for life to form and develop. Unfortunately, the magnetic 
measurements from the Cassini orbiter (Dougherty et al., 2004) 
do not reach the sensitivity required for a similar detection at 
Titan. The reasons for that are manifold, namely: (i) a fainter signal, 
since Saturn’s magnetic field is weaker at Titan’s orbit than that of 
Jupiter at its moons, (ii) a larger safe flyby distance of the 
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spacecraft above the surface, and more specifically, (iii) the shield- 
ing effect due to Titan’s ionosphere. 

Thanks to Titan’s dense atmosphere, the presence of a buried 
ocean might be detected by measuring the non-synchronous 
rotation of the crust (i.e., the length of day) due to the periodical 
exchange of angular momentum between the atmosphere and 
the crust (Tokano and Neubauer, 2005). An initial analysis of the 
location of geological features at different times led the Cassini 
Radar team (Stiles et al., 2008; Lorenz et al., 2008) to propose a 
non-synchronous rotation rate of 1.1 x 10~ 3 degree/day that 
would have been consistent with this model. However, the discov- 
ery of an error in the data reduction led to reconsider the rate as 
being rather 3 x 10 -4 degree/day (Stiles et al., 2010). Since this 
value is strongly correlated with the precession, it lies within the 
error bar of a non-synchronous rotation. Furthermore, a more elab- 
orate model of angular momentum exchange (Tokano et al., 201 1 ) 
predicts much smaller displacement, with amplitudes of polar 
motion from 100 to 1000 m for crustal thickness of 70-20 km, 
respectively. Given that the spatial resolution of the radar is several 
times larger than the 300 m surface sampling, it appears that the 
Cassini Radar observations could detect the non-synchronous com- 
ponent related to the presence of a buried ocean if the thickness of 
the crust is less than 20 km (O. Karatekin, Royal Observatory of 
Belgium, private communication). 

independently, the possibility to diagnose the shallow interior 
of planets, using the properties of Extremely Low Frequency 
(ELF) waves trapped within their atmospheric cavity, was sug- 
gested by Sentman (1990a) as a way to reveal the liquid metallic 
hydrogen surface of Jupiter. To date, the only observation that 
seems consistent with the presence of a subsurface water ocean 
at Titan is that of the ELF signal collected by the Huygens Probe 
(Beghin et al., 2010), as it descended through the atmosphere on 
January 14, 2005. 

The Permittivity, Wave and Altimetry (PWA) experiment (Grard 
et al., 1995), a subsystem of the Huygens Atmospheric Structure 
Instrument (HASI) (Fulchignoni et al„ 2002), was designed to per- 
form measurements of electron and ion conductivities, ELF and VLF 
(Very Low Frequency) waves during the descent through Titan’s 
atmosphere (Lebreton et al., 2005). One of the most exciting but 
speculative goals of the PWA experiment dealt with the possible 
existence of lightning activity in Titan’s atmosphere and the asso- 
ciated Schumann resonance (Schumann, 1952; Sentman, 1995), 
hereinafter referred to as “SR”. Subsequently, the presence of a 
ground or subterranean conductive surface, such a liquid water 
ocean predicted by theoretical models (e.g., Grasset and Sotin, 
1996), was expected to be a major scientific return of the mission. 
Preliminary data analyses revealed indeed that the 36 Hz signal 
observed with the PWA instrument (Grard et al., 2006) might sat- 
isfactorily fit the range of values predicted by models to be the sec- 
ond eigenmode of the SR (Nickolaenko et al., 2003; Simoes et al., 
2007). Given the absence of any acknowledged lightning activity 
on Titan after 72 flybys from 2005 to 2010 (Fischer and Gurnett, 
2011), it was proposed that the ELF turbulence generated within 
the ionospheric currents driven by the co-rotating plasma flow of 
Saturn’s magnetosphere (e.g., Sittler et al., 2009) might be the 
power source able to excite a “Schumann-like” resonance (Beghin 
et al., 2007). After reanalyzing the characteristics of the 36 Hz sig- 
nal in term of height profile and spin modulation, including a few 
attempts at identifying possible source regions from the measure- 
ments performed by the plasma and wave instruments on board 
the Cassini orbiter, the proposed powering mechanism was further 
developed in follow-up articles (Beghin et al., 2009, 2010). 

In the latter papers, because of the absence of any detectable 
radial electric field component of the 36 Hz signal, we were led 
to consider a mode of propagation with horizontal polarization. 
This assumption must be revisited and we consider instead the 


properties of waves with mixed horizontal and radial polarizations, 
called Transverse Electric and Magnetic (TEM) modes, known to be 
those of Earth’s SR (Wait, 1962). However, in a stratified cavity, 
bounded by an icy crust instead of a conductive ground and a mul- 
ti-layered profile of atmospheric conductivity, the TEM waves are 
known to degenerate into Longitudinal Section Magnetic (LSM) 
modes (Collin, 1991), whose properties are shown below to be 
notably different from those of Earth’s SR. 

We first establish the basic equations relevant to Titan 
conditions in Sections 2 and 3, referring to previous works and to 
Appendix A for some mathematical developments. After deriving 
the expressions of field components in Section 3 and the modal 
equations in Section 4, the formalism is subsequently applied to 
the PWA observations in Section 5 in order to constrain the param- 
eters of the cavity. Making use of the Huygens Probe dynamics data 
deduced from the Descent Imager Spectral Radiometer (DISR) 
observations (Karkoschka et al., 2007), we retrieve in Section 5.4 
a few case-studies of electric field strength profiles which compare 
favorably with Huygens measurements throughout the descent. 

Section 6 contains discussions regarding the two main results of 
our study, the first one being the constrained upper atmosphere 
electron conductivity profile above the Galactic Cosmic Ray 
(GCR) layer at around 60 km revealed by the PWA-MIP (Mutual 
Impedance Probe) and PWA-RP (Relaxation Probe) instruments 
(Hamelin et al., 2007; Lopez-Moreno et al., 2008). We invoke a 
probable relationship with thin aerosol layers and we consider a 
possible correlation with the zonal wind profile measured by the 
Doppler Wind Experiment (DWE) during the descent (Bird et al., 
2005). The second major result concerns the crust characteristics 
and the presence of a buried conductive ocean, whose constrained 
parameters are discussed in terms of Titan’s interior thermal 
dynamics. Although only few ground truth and remote sensing 
measurements support the relevance of the constraints derived 
from our work, we are able to conclude in Section 7 that our inter- 
pretation of PWA-ELF data is consistent with recent models of Ti- 
tan’s interior. The cavity appears to be bounded upwards by thin 
ionized layers of aerosols, at least up to an altitude of 140 km 
(Lopez-Moreno et al., 2008) and downwards by a water-ammonia 
ocean buried beneath a semi-rigid icy crust, most likely 55-80 km 
thick. 

We revisit, in Appendix B, a few experimental open points that 
have been subject of recurring questions regarding the validity of 
the PWA measurements, such as: (i) the absence of multi-SR 
modes that should have been observed in addition to the 36 Hz 
line, (ii) the disappearance of the signal 16 s after landing, and 
(iii) a presumably artifact induced by boom-antenna vibrations at 
36 Hz possibly triggered by wind-blasts or by the descent velocity 
airflow (Beghin et al., 2007). 


2. Formulation of the problem 

Once compendious analyses of the raw data had failed to reveal 
any convincing signature of a radial electric field component 
(Beghin et al., 2009), it was wrongly concluded that the polariza- 
tion might be essentially horizontal. Such an erroneous conclusion 
came from the inaccuracy in the time resolution of the 36 Hz signal 
strength which did not exhibit convincing correlation with the fast 
variations of the tilt of Huygens vertical axis (Section 6.1), as could 
be expected from the presence of a significant radial electric field. 
When the excitation comes from a horizontal electric dipole, or the 
system of azimuthal current sheets proposed by Beghin et al. 
(2007), horizontal and radial polarized waves coexist (Bahar and 
Fitzwater, 1983), and trapped TEM modes have both radial and 
horizontal electric components, but no radial magnetic field com- 
ponent (Galejs, 1972). However, as the structure of Titan’s cavity 
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is stratified, we must instead consider the LSM modes which are 
known to exhibit high quality factors in wave-guides filled with 
dielectric layers (Collin, 1991). 

Though the following analytic treatment presents physical anal- 
ogies with that of TEM modes of Earth’s SR (Greifinger and Greifin- 
ger, 1978; Sentman, 1990b; Fiillekrug, 2000; Nickolaenko et al., 
2003), both approaches differ significantly, not only because iono- 
spheric horizontal current sources are involved instead of vertical 
currents, as would be the case for ground-to-cloud return-strokes, 
but also because of the presence of stratifications which requires 
appropriate treatments at the layer interfaces (Wait, 1962; Collin, 
1991). The physical concepts are nevertheless of the same nature, 
in the sense that we treat the scattering mechanism of an incident 
plane wave over a perfectly conductive inner surface, which is 
found to be a satisfying approximation compared to results of 
full-wave expansion treatments (Bahar and Fitzwater, 1983). The 
conductivity of Titan’s icy-crust (~0.2nS nr 1 ) measured by the 
PWA-MI instrument (Grard et al., 2006) is too low to reflect ELF 
waves, and its permittivity (~2.5) is found consistent with the va- 
lue derived from the Cassini Radar Scatterometer (Wye et al., 
2007). Therefore, we consider a subterranean conductive surface, 
presumed to be a water-ammonia liquid ocean buried under an 
icy crust made of dielectric material (Fig. 1 ). For the electrical prop- 
erties of the atmosphere, we consider the electron conductivity 
profile deduced from the PWA-MIP-RP instruments up to about 
90 km, including the GCR layer. At higher altitudes, the perfor- 
mance of PWA-MIP instrument is doubtful (Hamelin et al., 2007; 
Beghin et al., 2009), due to a suspected non-nominal configuration 
involving probably the boom-antenna release system prior to the 
main parachute jettison (Lebreton et al., 2005). Moreover, above 
the GCR layer, the PWA-RP measurements exhibited discontinu- 
ities attributed to stratifications of few hundred-meter size with 
the occurrence of “plateaus”, i.e., time intervals during which the 
response of the probe potential remains constant for 6-32 s, at 
about 100 and 72 km respectively (Lopez-Moreno et al., 2008). 
Thus, in the absence of high altitude measurements applicable to 
a global planetary scale, we will resort to case-study models with 
one-scale-height atmospheric electron conductivity yielding a best 
fit between retrieved and experimental field profiles (Sections 4 
and 5.5). 
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Fig. 1. Model of Titan’s layered cavity. Dash-dot line : altitude profile of |n 2 | 
measured by the PWA-MIP-RP instruments from 0 to 90 km altitude, and analytic 
continuation upwards ( dot line). Vertical scales: z is referring to inner boundary 
(ocean), h is the altitude above Titan surface; hi (zi) and h 2 ( z 2 ) are the conduction 
and diffusion altitudes, respectively (see text in Section 4). Values of each 
parameter are given in Table 2 for six case-studies. 
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In a spherical geometry the radial structure is assumed to be ar- 
ranged in series of discrete shell slabs wherein the electrical prop- 
erties are homogeneous, so that we may apply the following usual 
Maxwell equations 

VxH=J ( +® V ■ B = 0; VxE = -^ ; D = e 0 n 2 E: B = 

( 1 ) 

where H and E are the magnetic and electric vectors, D and B the 
electric displacement and the magnetic induction; e 0 and /r 0 are 
the free space permittivity and permeability, and J e is a possible in- 
duced surface-current distribution, as for instance, at the inner per- 
fectly conductive boundary. The complex relative permittivity e of 
the medium is written in terms of conductivity as 

£ = n 2 = Ree + i Im s, with 1m e = , (2) 

CO £ 0 ' 

where n is the index of refraction and a the scalar conductivity. All 
quantities are assumed to vary versus time as exp (- icot ), where /is 
the wave frequency, to = 27t/and i = ^/-l. Introducing the Riemann- 
Lorenz vector potential A, from the second and third equations in 
Eq. (1) one readily obtains 


B = V xA 


and V x 



(3) 


Introducing the scalar potential function V F, the second equation in 
Eq. (3) yields 

f) a 

V‘f' + E + ^ = 0. (4) 

at 

Applying the curl operator to the first equation in Eq. (3), and intro- 
ducing the identity k 2 = <±> 2 [i 0 £ 0 , where k is the free space wave 
number defined as 2 n/X, and X is the free space wavelength, we 
get the wave equation as a function of both potentials 


], 2 n 2 

V x VxA=i— Vf + lfVA. 
m 


(5) 


Using a procedure similar to that applied for the computation of ter- 
restrial SR parameters (Sentman, 1990b), the wave Eq. (5) may be 
reduced into a set of two independent differential equations for 
each potential A and *P. The procedure, developed in Appendix A, 
is basically the same as for radially polarized TM modes, i.e., when 
the radial magnetic field is forced to zero, with the consequence for 
the vector potential to have only a radial component A r . The result- 
ing Eqs. (A8a) and (A8b) are known as particular forms of general- 
ized Helmholtz wave equations applied to inhomogeneous cavities, 
and will be referred to in the next sections whenever needed. 


3. Field expressions and approximate solution of wave equation 

First, we define the initial boundary conditions for the horizon- 
tal components of electric and magnetic fields obeying Eqs. (3) and 
(4), i.e., 


{E„ = E„ = 0} z=0 > [0 VI', 0: <p V'f'r = 0; r ■ V^ r #0; v F r = 0} z=0 , 

(6a) 

{H,)#0; H^O } z=0 => {A r = A 0 } z=0 , (6b) 

where A 0 is the amplitude of the vector potential at the lower con- 
ductive boundary (z = 0) after performing the change of variable 
z = r - a + z c , where a is the mean Titan’s radius and z c the thickness 
of the crust (Fig. 1 ). Starting from the wave equations Eqs. (A8a) and 
(A8b), we may proceed either with the vector potential or the scalar 
potential, independently of each other. The field expressions 
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deduced from Eqs. (3), (4), (A3), and (A7), are written in terms of the 
vector potential as 


H,i = 0- 


A r dY(8, cp) 


A r dY(9, cp) 


a)l a sin 0 dcp 
"io>l(l + 1) /A 

E p = (p- 


k 2 a 2 
ico 


Y{B,cpy, e ( ) = e 


a)i 0 dd 


k a V 

dA r ~ 

Oz 


1 0A r \ dY 
k 2 a sin ff{n 2 dz ) dcp 


(7a) 


(7b) 


Note that the LSM modes are entirely defined by a single scalar 
function A r , namely the Hertzian potential (Collin, 1991). Similarly, 
the fields can be written in terms of the scalar potential 


(ArY)z- 


(ArY) z 


\£ c dz ) z - 


J_aV\ 

n l dz ) z : 


(13) 


where zy and z,f refer to the levels below and above the crust sur- 
face, respectively. Applying to Eq. (11) the continuity conditions 
for A r and its first derivative, we get the wave expression that is va- 
lid upwards throughout the atmosphere as long as each slab is sup- 
posed to be homogeneous 


(Ar) z> z c — 


COS(ffiZc) cos I< 2 (z - z c ) 


njK i 

e c I < 2 


sin(K,z c ) sin I< 2 (z 



(14) 


H e -ik 2 n 2 a d'V, dY 

[k 2 n 2 a 2 -1(1+ l)]/( 0 sin0 dz dcp 
ik 2 n 2 a d^r dY 
tp ~ <P [k 2 n 2 a 2 - 1(1 + 1 )}ju 0 ~dz dd ’ 


(8a) 


E r = r 


1(1+1) 


0% 


k 2 n 2 a 2 - 1(1 + 1) dz 


Y(9,cp ); E„ = 0 


-y r dY 

~a~d8' 


— cp 


-¥ r dY 
asin0 dcp' 


(8b) 


We may derive an approximate solution of Eq. (A8b) for the Hertz- 
ian potential, starting from the lower boundary z = 0 up to the ion- 
ospheric boundary at the altitude fq (Fig. 1). The latter is defined as 
the conduction boundary (Greifinger and Greifinger, 1978), where 
the conduction current becomes equal to the displacement current. 
According to Eq. (2), this condition reads 

<= 1+ W Sl+i - (9) 


In the atmosphere up to /q, the electron conductivity measured by 
PWA-MIP-RP instruments is sufficiently low to consider |n 2 | close to 
1. Applying the wave equation within slabs, where n 2 is assumed to 
be constant, Eq. (A8b) simplifies as 


W9AA 

n 2 \dzj 


~ A) 


^cos(/C,z c )sin/<: 2 (Z-Z c ) 

- n 2 


+ — sin(/C 1 z c )cos/C 2 (Z-Z c ) , 


(15) 


where K 2 is the complex radial wave number solution of Eq. (10) in 
the atmosphere when n\ ss 1. For / = 2, /= 36 Hz and a ~ 2600 km it 
can be seen from Eq. (10) that y 2 ~ 1.56; thus I< 2 is almost purely 
real negative as long as n\ ~ 1, so that K 2 is lying nearly on positive 
imaginary axis. This fact does not imply that the wave is damped, 
because the phase increment \K 2 z\ through the atmospheric portion 
of the cavity is nearly zero. Thus, the first order expansions of Eqs. 
(14) and (15) are real, so that A r does not depart significantly from 
its initial value A 0 at least up to the ionospheric boundary hi, 
although its first derivative 0A r /Sz is not zero. As a consequence, 
the horizontal magnetic field components keep nearly constant 
amplitude versus z for given 0 and cp by virtue of Eq. (7a). Introduc- 
ing Eq. (15) in the last two expressions of Eq. (7b) one can check 
that the horizontal electric components increase linearly from 0 at 
z = 0, up to their maximum value reached in the vicinity of fq, while 
apart from the GCR layer, the radial electric component is almost 
constant from the ground up to around /q (Section 5.3). 


4. Modal equation 


° Dz 2 + K 2 A r = 0 with K) = k 2 (nj - y,j and y, = , (10) 

where j is the slab index, starting from j = 1 at the bottom of the icy 
crust (z = 0), Kj is the radial wave number, / is the SR eigenmode de- 
gree. As long as we consider ELF waves with free space wavelengths 
X of the order or even larger than a, the value of Kj in the Eq. (10) is 
of the order of magnitude of (27t/A) 2 , thus |/( 2 z 2 | <c 1 everywhere in 
the cavity for Kj positive, negative or complex. Although the general 
solution of Eq. (10) is known to be a combination of Bessel functions 
for the complex variable Kj, we consider here the exponential thin- 
shell approximation (Galejs, 1972), whose solution is a combination 
of upgoing and downgoing harmonic waves, as 

A r (z) = a exp(iKjZ) + fl exp(-iKjZ), (11) 

where a and fl are determined by the initial boundary conditions 
Eqs. (6a) and (6b), i.e„ a = I) = A 0 I2. Thus, the solution of Eq. (11) 
within the crust yields 

A r (z)=A 0 cos(ff 1 z); i^ = -A 0 ^sin(KiZ) forz<z c , (12) 

c>c C/Z oq 

where K, is given by Eq. (10) with n 2 = e c . We apply now the inter- 
face conditions to each crossing of successive slabs, in such a way 
that the radial component of the displacement vector D and both 
horizontal electric and magnetic field components are continuous. 
Starting from the ground surface, this reads 


The modal expression of the eigenmodes of Titan’s cavity is de- 
rived using the approach outlined by Greifinger and Greifinger 
(1978) for the terrestrial SR in cylindrical geometry and later on 
in spherical geometry by Sentman (1990b) and for Titan by 
Nickolaenko et al. (2003). Here, the only significant difference with 
respect to the Earth conditions is the presence of the subsurface 
icy-crust, whose physical parameters must be included in the mod- 
el. We saw that the amplitude of A,-, given by Eq. (14), is almost 
constant as long as |/C 2 z| < 1. Thus, we may perform the integra- 
tion of the height-gain Eq. (A7) for ¥ r throughout the cavity up 
to some distance above hi, where A r may be taken outside the 
integral. Inside the icy crust, the result is straightforward 

Vr (z) = icoA 0 (1 - y,/s c )z. (16) 

Boundary conditions similar to Eq. (13) applied to the scalar poten- 
tial read 


V r j~V rJ _i + (z 



and 


(rfdVr) 

[Kj dz ). 


n U 

d'Vr 




for Zj ! < z < Zj, 


(17) 


where j and j - 1 are the indices respectively above and below any 
slab boundary. It can be checked from Eq. (17) that *F r is continuous 
across the ground surface. This is not the case for its derivative 
which exhibits an abrupt 90° phase shift while the permittivity 
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drops from e c > 1 down to 1. However, according to Eqs. (8b) and 
(17), the horizontal electric components grow linearly in the atmo- 
sphere, as we found above using the vector potential (Section 3). For 
simplification purposes, we choose the analytical approximation of 
one-scale-height exponential profile which is shown to predict sat- 
isfactorily the Earth’s SR modes for normal ionospheric conditions 
(Fullekrug, 2000). The presence of small irregularities, such as the 
GCR layer (Fig. 1), could introduce minor deviations which are 
found to have a weak effect on field profiles (Section 5.3). The 
one-scale-height exponential profile referring to the altitude of con- 
duction boundary h t , where n 2 = 1 + i, reads 

n 2 = 1 + i exp , (18) 

where f is the scale height and z^ =z c + h ^ (Fig. 1). Using Eq. (18), 
and integrating Eq. (A7) from the ground up to z above Zi, yields 

*Pr(Z) = 'Pr(Zc) + icoAo jf (l - J)dz 

= ¥ r (z c ) + icwAo[z c (y,-l)+z-y,(zi (19) 

Substituting in Eq. (19) the value of 4V (z c ) given by Eq. (16), we 
obtain 


¥r{z) = icoA 0 [z-y l (z l -Zc+z c /Sc-i^)]; (z>Zi+f). (20) 

In regions well above z lt n 2 reaches high values such as the last term 
with gradient of rr 2 in Eq. (A8a) vanishes. Similarly, we may ignore 
the term y t compared to n 2 in the coefficient of 4V. Then, using 
Eq. (9) and the identity k 2 = /( 0 e 0 co 2 , the wave equation reduces to 

^y^ + V( 0 w<r(h)4V = 0; (z>Zi + 0. (21) 

We define the upper ionospheric boundary as z 2 = h 2 + z c , where h 2 
(Fig. 1 ) is the altitude of diffusion for the magnetic components and 
the ultimate reflection boundary for the horizontal electric compo- 
nents (Greifinger and Greifinger, 1978). Extrapolation of the con- 
ductivity profile given by Eq. (18), yields the expression of the 
conductivity at z = z 2 

er(z) = a 2 (z 2 ) exp with er 2 (z 2 ) = e 0 coexp 

( 22 ) 


z 2 -z, 

C 


Performing the change of variable v = exp[in/4 + (z - z 2 )/20, Eq. 
(21) becomes 


d 2 4V 

d 2 v 


+ 


1 d'V r 

V dv 


+ 4/i 0 a>er 2 ( 2 4V = 0. 


(23) 


Giving the conductivity <r 2 a value such that 4)i 0 mcr 2 £ 2 = 1 at the 
boundary h 2 , Eq. (23) takes the form of a Bessel’s differential equa- 
tion whose solution is proportional to the zeroth order Hankel’s 
function of first kind H' n ' ' (v) (Abramowitz and Stegun, 1972). Above 
z 2 the function tends rapidly towards zero when u -*■ oo exp in/4 for 
z » z 2 , whereas below z 2 , where i)-*0 exp in/4 for upgoing waves, 
the small argument expansion of Hankel’s function takes the form 


relation 4f.i 0 u>a 2 i f = 1, we may link both boundaries through the 
scale height, as follows 


h 2 


h, +£ln— y-. 
4 k 2 C 2 


(26) 


As a consequence, Eq. (25) involves only five unknown parameters 
of the cavity, namely: (i) the thickness of the icy crust z c , (ii) the real 
part of the crust permittivity Re e c , (iii) its loss tangent <5, assumed to 
be weak enough such as S -c Re e c , (iv) the altitude of the conduction 
boundary h lf and (v) the scale height of the conductivity 0 Consid- 
ering f <tc h, and merging Eq. (25) with the expression of yi in Eq. 
(10), it is straightforward to derive the approximate modal equa- 
tions in terms of eigenfrequency and quality factor as 


coi = 


1(1 + V 


h i + Zc/Ree c 


1 1/2 


1 — i 


hi +z c + fln(l/4/< 2 C 2 ) 

. ( z c < VRee c + 7tC/4 ntj/4 

V ^i+z c /Ree c + /ii+ Zc -Kln(l/4/< 2 0) 


(27) 


Q = 


Re co 
2 Im co 


2z c <5/Ree c + 710/2 710/2 

hi+z c /Ree c ^ +z c + £ln(l/4fc 2 £ 2 ) 


(28) 


One may check that for z c = 0 the above expressions reduce to the 
usual ones for the terrestrial SR (Fullekrug, 2000). It is worth 
emphasizing the crucial contribution of crustal parameters for con- 
trolling the Titan’s SR eigenmodes, which makes a fundamental dif- 
ference with the Earth’s SR. In the next sections we take advantage 
of this particularity for constraining the cavity parameters. 


5. PWA observations and constrained model cavity 

The average profiles of atmospheric electron conductivity and 
36 Hz signal strength measured with the PWA-MIP-RP-ELF instru- 
ments (Hamelin et al., 2007; Lopez-Moreno et al., 2008; Beghin 
et al., 2009) are plotted in Fig. 2, in middle-left and middle-right 
panels, respectively. For the purpose of comparison we have in- 
cluded a few additional measurements. The first one (left panel) 
is the zonal wind speed profile obtained with the Doppler Wind 
Experiment (DWE), conducted with the Green Bank and Parkes 
Radio Telescopes (Bird et al., 2005), plotted along with the Huygens 
Probe descent velocity. The last right panel plots are the raw and 
smoothed tilt variations of the Huygens vertical axis derived from 
the DISR experiment (Karkoschka et al., 2007). Prior to proceeding 
with the analysis of PWA data we address the efficiency of the ELF 
emission mechanism for further comparison with the strength of 
the Huygens signal. 

5.3. Modulation mechanism of the current source 


4 'r(h) = c(z-z 2 -i^c) for z < z 2 , 


(24) 


where C is an arbitrary coefficient. Matching Eqs. (20) and (24) 
within the common region z x < z < z 2 implies that C = icoA 0 , then 
the modal equation reads 


Z 2 + 1071/2 

Zi — z c (l — l/e c ) — 1071/2' 


(25) 


It becomes convenient here to use the altitude coordinate h = z- z c 
(Fig. 1) while referring to atmospheric profiles. Using the above 


The current-driven ion-acoustic instability mechanism pre- 
sumed to be the ionospheric source of Titan’s SR (Beghin et al., 
2007), also called two-stream instability, was known long ago to 
develop in space and laboratory collisional plasmas (e.g. Bychen- 
kov et al., 1994). The instability is due to the anomalous conductiv- 
ity that occurs as soon as the relative drift velocity between 
electrons and positive ions exceeds the local sound speed, and 
the ELF modulation bandwidth stretches up to the ion-plasma fre- 
quency (Appendix B.l). Since the modulated currents are draped 
around Titan’s ram hemisphere (Ledvina and Cravens, 1998), they 
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Fig. 2. Smoothed altitude profiles of data used in this study. From left to right, (i) zonal wind-speed ( solid line), from DWE (Bird et al„ 2005), and descent velocity (dot line), 
from the mission event list; (ii) average electron conductivity from PWA-MIP-RP; (iii) smoothed field strength of the 36 Hz line from PWA-ELF; (iv) raw and smoothed profiles 
of the tilt of Huygens vertical axis, from D1SR (Karkoschka et al„ 2007). Note the coincidence between the wind shear above 60 km and the GCR layer (horizontal dot line). 


may be compared to vibrating strings sustaining the resonance of 
the cavity. 

On Earth, the natural excitation of the ion-acoustic instability 
has been observed locally in the equatorial electrojet by rocket 
borne instruments in form of electrostatic wave at 25 Hz with 
2 mV itT 1 amplitude (Pfaff et al., 1987). Likewise, the current mod- 
ulation mechanism was proved to be an efficient stimulus for ELF 
waves with an active experiment carried out in the region of the 
terrestrial polar electrojet (Me Carrick et al., 1990), using the 
High-Power Auroral Stimulation (HIPAS) HF (High Frequency) 
heating facility, near Fairbanks, Alaska. In this experiment we as- 
sume that the generation of ELF waves was coherently stimulated 
by forcing the anomalous conductivity of the medium with depo- 
sition of modulated powerful HF waves, in presence of the polar 
electrojet. Electromagnetic ELF waves (6-76 Hz) with magnetic 
field amplitudes of up to 2 pT and 0.2 mV m _1 electric field were 
received at the ground, 35 km from the HIPAS site. Considering 
that the E region illuminated with the HF source was covered by 
a large number of sidelobes distributed uniformly in azimuth over 
45° (Me Carrick et al., 1990), the extent of the disturbed region was 
roughly 90 km, to be compared to about 60° of latitudinal extent 
(~2800 km) for Titan’s current sheets. Applying the geometrical 
scale factor 90:2800, one deduces the transfer ratio 2:62 pT that 
will be considered further in Section 5.4. Notice that this ratio is 
probably underestimated, because the HIPAS fixed frequency stim- 
ulus was coherent over 90 km, contrary to random noise of Titan’s 
ion-acoustic turbulence. 

5.2. Mode degree and subsurface parameters constrained by PWA 
observations 

If we assume that Titan is an ideal Schumann cavity with a 
nominal radius of 2575 km, introducing yi = 1 in Eq. (10), the first 
three eigenfrequencies are found roughly 26, 45, 64 Hz for 1 = 1, 
2, 3 respectively. The frequencies predicted by models are either 
close to these values (Nickolaenko et al., 2003) or widely spread 
around, depending on the cavity structure (Simoes et al., 2007). 
Thus, the first parameter to be confirmed is the degree / of the 
36 Hz signal. Assuming that the coefficient is constant for all 
modes, introducing /= 36 Hz, / = 2 in Eq. (10) with a = 2600 km 
including partly the dielectric crust contribution, we get y ; = 1.56. 
From Eq. (25) it appears that y ( cannot be lower than unity since 
z 2 > Zi by definition (Fig. 1 ). Thus 36 Hz cannot be the fundamental 
eigenmode. However, this frequency could be either the second or 
the third mode (Simoes et al., 2007). Thus, we need additional data 
to determine unambiguously the mode degree. It is a common fea- 
ture of exponential scale-height ionospheric models (e.g., Sentman, 
1990b), that both horizontal and radial electric components reach 


maximum amplitudes at around the conduction boundary One 
can see in Fig. 2 that the maximum strength of the 36 Hz signal lies 
just below 100 km altitude. Therefore, whatever combination of 
electric field components might have been measured by PWA an- 
tenna (Section 5.3), we may constrain lq to lie at around 100 km. 

The third constraining parameter is the quality factor Q~6 
(Fig. 3, left-hand panel) derived from an average of 61 individual 
ELF power spectra with 16 frequency bins each, measured between 
the altitudes of 90.7-80.5 km, where the 36 Hz signal strength is 
maximum. Although this value is comparable to the terrestrial 
ones (e.g., Chand et al., 2009), it is almost one order of magnitude 
larger than predicted by Titan’s models with a perfectly conductive 
ground (Nickolaenko et al., 2003), whereas Q values as large as 5 
are possible in the presence of a perfectly conductive surface be- 
neath a crust 100 km thick (Simoes et al., 2007). Thus, introducing 
Q_= 6 in the modal Eq. (28) allows us to cross-constrain the remain- 
ing four unknown parameters, namely: the thickness z c of the 
crust, its complex permittivity (Re e c and loss tangent <5), and the 
conductivity scale height f of the upper boundary. The loci the ma- 
trix solutions are plotted in Fig. 4 for the two options / = 2 and 3, as 
a function of the scale-height (horizontal axis) and the crust per- 
mittivity (vertical axis) for different values of crust thickness. 

The only available measurement of the surface permittivity in 
the ELF range comes from the PWA-MIP data at the Huygens land- 
ing site. The sounding area embodied about one meter depth and a 
first estimate led to Re e c ~ 1.8, subsequently constrained between 
1.9 and 3.2 (Grard et al., 2006). In the microwave range of the Cas- 
sini Radar Scatterometer, the mean ground permittivity of Titan’s 
surface is Re e c = 2.2 ± 0.05 (Wye et al., 2007). However, it is not 
obvious that the electric characteristics of the crust should be the 
same at depth, because the ice seems to be contaminated by 
deposits of organic residue such as tholins (Zarnecki et al., 2005). 
In a crust made of water Ice I (Tobie et al., 2005) with chemical 
inclusions such as ammonia concentration up to 10%, the permit- 
tivity might reach 3-4 (Lorenz, 1998), and may vary with increas- 
ing temperature at depth (Section 6.2). Assuming an unlikely wide 
range of crust permittivity 2 < Re e c < 5, and ionospheric conductiv- 
ity scale height 3 < £ < 15 km, for / = 3 we find that the thickness of 
the crust should be 400-800 km (Fig. 4, left panel). This is much 
larger than the upper bound (~160 km) predicted by recent mod- 
els (Deschamps et al., 2010). Moreover, even for models with unre- 
alistic low amount of radioactive elements in the silicate fraction 
(e.g., Grasset and Sotin, 1996), the presence of ammonia within 
the crust is not compatible with such thickness because a low eu- 
tectic temperature leads to partial melt and formation of an ocean 
at a maximum depth of about 100 km. Thus, the mode 1 = 3 is 
incompatible with a frequency of 36 Hz and a quality factor of 6. 
On the contrary, for / = 2, a permittivity 2 < Re e c <4 and loss 
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Fig. 3. Asterisks, left panel : average of 61 PWA-ELF power spectra data in 3 Hz resolution mode, without the missing Channel A odd bins [(2n -1)3]. Right panel : power 
spectrum (even bins [(2n) 6] Hz) transmitted 8 s after touchdown in 6 Hz resolution mode. Dot lines : smoothed background noise. Shaded areas: power spectra expected to 
have been observed with nominal frequency resolution. Vertical solid lines in left panel: quadratic mean composite of expected side band values plus noise fitting the 
measurements ( asterisks ). 
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Fig. 4. Loci of cross-constrained crustal parameters versus the ionospheric scale-height for eigenmodes / = 3 ( left panel) and l = 2 ( right panel). 


tangent 0.25 < S <0.5 (i.e., a conductivity a =1 to 4nS m _1 ), the 
crust might be 40-80 km thick, and the ionospheric conductivity 
scale-height between 3 and 4 km (Fig. 4, right panel). Thus, the 
36 Hz frequency fulfils the required conditions to be the second 
eigenfrequency of Titan’s SR, denoted below as an LSM 2 mode. 

5.3. Altitude profiles of field components of the LSM 2 mode 

Using the formalism developed in Sections 3 and 4, our purpose 
was to check whether the constrained parameters of the cavity are 
compatible with the height profile of the PWA signal strength. The 
first step of the procedure consisted in computing the field profiles 
for a few case-studies considering the most likely values of the 
constrained parameters. The second step consisted in retrieving 
the profiles of the electric signals that would have been collected 
by the PWA antenna, and to compare afterwards the retrieved pro- 
files with the measurements. The model field profiles are com- 
puted integrating the scalar potential Eq. (A7) through successive 
slabs 5 km thick along the vertical path outlined in Fig. f. At each 
interface we apply the transfer conditions Eq. (17) taking care of 
their analytic continuity within the complex plane, namely, the 
sudden phase shift at the ground-atmosphere interface mentioned 
in Section 4. The normalization factor (coA 0 ) is given an arbitrary 
value of 1 Virr 1 , while A r is kept constant up to the vicinity of z 2 . 


The field amplitudes are derived from Eqs. (8a) and (8b) using 
the usual spherical harmonic function 

Y(8, tp) = P["(c os 9) cos mtp, (29) 

and its derivatives, where P[" are Legendre polynomials of degree / 
and order m, the latter being an integer or zero. The terrestrial ver- 
tical current sources are suited to traditional SR modes which are 
characterized by a single azimuthal magnetic component and 
m = 0 (Sentman, 1990b; Nickolaenko and Sentman, 2007). On the 
contrary, Titan’s horizontal current sheets (Fig. 5) constrain m to 
be an integer (Wait, 1962). We must, therefore, consider the tesser- 
al modes LSM[" to denote Titan’s SR, with m = 1 for the fundamental 
order. For l = 2, the Legendre polynomial of first order reads 
P 2 = -3/2 sin 2 0. Introducing this expression in Eq. (29) we obtain 

3 

Y(8, cp) = -- sin 28 cos tp, 

d Y d Y 3 ( 30 ) 

— = -3 cos 20 cos tp ; — = - sin 28 sin tp, 
a0 dtp 2 

where 0 and tp are the colatitude and the longitude of the observer 
with respect to the source. 

Without loss of generality, we consider the configuration 
plotted in Fig. 5, where the source and the observer lie within a 
reference system with an ideal north-south axis orientation of 
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Fig. 5. Spherical coordinates for arbitrary observer’s location (M) with respect to a 
horizontal elementary current dipole Ids, in the ideal configuration of azimuthal 
current sheet. 


the co-rotating Saturn’s dipolar magnetospheric plasma flow (Led- 
vina and Cravens, 1998). For better understanding of the differ- 
ences between Titan’s horizontal and Earth’s vertical current 
sources, the expressions of Y(0, tp) for the modes 1 = 1, 2, and 
m = 1,0 respectively, are reported in Table 1. It is seen that the ant- 
inodes (peak values) of the radial electric component for the sec- 
ond Titan’s tesseral mode LSMj take place at latitudes ±45° 
(±180°) with respect to the source region, contrary to the Earth’s 
TEIVl” zonal mode, where the antinodes lie at 0° and 180°. 

Introducing Eq. (30) in the field expressions Eqs. (8a) and (8b), 
the numerical integration of Eq. (17) was performed by steps of 
5 km starting from z = 0 up to z 2 (Fig. 1), for constrained parame- 
ters derived from Fig. 4 (right panel) and different values of 0, cp. 
A few case-studies for z c ranging between 40 and 80 km are re- 
ported in Table 2, and the retrieved field profiles for case-study 
#3 are plotted in Fig. 6. With the exception of the weak contribu- 
tion of the GCR layer, the altitude profiles of the horizontal compo- 
nents are similar to those of the terrestrial two-scale-height model 
(Sentman, 1990b). The slope of the horizontal electric component 
(|£ h |) profile does not change substantially across the ground sur- 
face, reflecting the fact that the conductivity of the crust is low 
in this example (tan 5 = 0.27, i.e., a = 1.34 nS nr 1 ). However, the 
profile of the radial electric field component exhibits a different 
behavior with respect to Earth’s conditions. The amplitude of £ r 
rises suddenly across the ground surface (Fig. 6) while the permit- 
tivity drops from e c to 1, as it can be checked by substituting the 
interface relationship for 8q>/Sz given by Eq. (17) in the expression 
of £ r in Eq. (8b). 


5.4. Retrieved amplitude profile of the measured electric field 

The next step consisted in retrieving the Huygens measure- 
ments obtained with the above model field profiles. The PWA 


electric antenna was made of two short monopoles, each 45 cm 
long, mounted on two opposite sides of the gondola (Grard et al., 
2006), yielding an effective length of about 1.6 m along the Huy- 
gens horizontal axis (Cadene, 1995). Such conditions imposed by 
severe constraints of the mission are obviously not comparable 
with several meters long vertical antennas installed onboard 
Earth’s balloon flights (Ogawa et al., 1979). However, since the re- 
sponse to electromagnetic fields of a spinning horizontal antenna 
is attitude dependent, the modulation of the 36 Hz signal has been 
revealed by a statistical analysis of time intervals comparable to 
the spin period (15-35 s) of the Huygens gondola (Beghin et al., 
2009). In the present study, in order to reveal the global contribu- 
tion of all electric components associated with the variations of the 
tilt of Huygens vertical axis, we have computed the signal strength 
profiles that should have been measured in the presence of model 
field profiles such as those plotted in Fig. 6. 

The measured signal strength is the root mean square of hori- 
zontal E h and radial £ r components projected along the antenna, 
which reads 

£m = [(E/i cos Q) 2 + (£ r sin f?) 2 ] 1 ^ 2 with Q = [1 - cos 2 i/cos 2 t] 1/2 , 

(31) 

where Q is the smoothed tilt value (Fig. 2, right panel), computed 
from values of the roll (p) and pitch (t) angles (Table 1, in Kar- 
koschka et al., 2007). The retrieved profiles, corresponding to three 
case-studies (#2, 3, 4) from Table 2, are plotted along with the 
experimental profile in Fig. 7. While the altitude hi of conductivity 
threshold has been prescribed at 100 km in the model, the peak val- 
ues of retrieved responses are observed at 95 ± 5 km (i.e., one step 
of the model), fitting closely the spread of experimental values. This 
satisfactory model-data fit allows us to confirm the altitude of hi at 
around 100 km for the model of atmospheric conductivity discussed 
in Section 5.5. 

Note that the retrieved profiles in Fig. 7 have been plotted using 
a homothetic transformation scale in order to match their peak 
amplitude with the experimental value, so that the scale ratio 
yields a rough estimate of the absolute field amplitudes plotted 
in Fig. 6. It appears from Eq. (31) that the main contribution to 
the 36 Hz signal measured at around 90 km comes from the dom- 
inant component £ r while the average tilt plotted in Fig. 2 (right 
panel) at this altitude is about 9°. Substituting the latter value 
and the measured amplitude |£ M | = 5.5 mV m 1 Hz~ 1/2 in Eq. (31), 
neglecting the amplitude of E h which is 40 times smaller than £ r 
at this altitude (Fig. 6), we obtain an estimate of the power spectral 
density of the radial electric component |£ r | ~ 35 mV m _1 Hz~ 1/2 at 
90 km. This value must be compared with the arbitrary amplitude 
of 1600 mV nr 1 (Fig. 6) obtained with the normalization factor 
ojA 0 = 1 Vnr 1 . According to the experimental calibration ratio 
between sine amplitude and power spectral density, applying the 
homothetic ratio (35:1600) to all components plotted in Fig. 6, 
the amplitude of the horizontal electric component £,,| is found 
~ 1 mV nr 1 at around h i, and the horizontal magnetic induction 


Table 1 

Top three rows: Expressions of the spherical harmonic function Y(9, cp) and its derivatives for the first two SR modes of Titan and Earth: (i) Titan’s tesseral modes for horizontal 
current source, (ii) Earth’s zonal modes for radial current source. Bottom two rows: Latitudes of antinodes 0 M (|Er|Max) and nodes 9 m (|E r | m i n ) of radial electric field components with 
respect to the source location. 



Titan’s tesseral modes 


Traditional Earth’s zonal modes 


LSMj(/= l,m = 1) 

LSM)(( = 2,m = 1) 

TEM?(/ = l,m = 0) 

TEM^(/ = 2,m = 0) 

Y(9, cp) 

-sin 9 cos cp 

-3 sin 9 cos 6 cos (p 

cos 9 

Yt (3 cos 2 9- 1) 

dY(6, (p)ldB 

-cos 6 cos (p 

-3 cos 29 cos (p 

-sin 9 

-3 sin 9 cos 9 

dY(B, (p)ldtp 

sin 9 sin cp 

3 sin 9 cos 9 sin cp 

0 

0 

(|E r |lVlax) 

±90° 

±45° (±180°) 

0°, 180° 

0°, 180° 

Bm (|Er|min) 

0°, 180° 

0°, 180° 

±90° 

±54.7° (±180°) 
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Table 2 

Case-studies of Titan's cavity parameters constrained by the plots in Fig. 4 (right panel). 


Parameters 

z c (km) 

Re £ c 

tan 8 

hi (km) 

Co -60 (km) 

^60-80 (km) 

( 80 -h 2 (km) 

h 2 (km) 

a 2 (Sm 1 ) 

0 (deg) 

9 (deg) 

Case-study 

#1 

40 

3 

0.37 

100 

+5.2 

-5.2 

3.5 

136 

7 x ur 5 

50 

20 

#2 

50 

2.65 

0.295 

105 

+5.2 

-5.2 

3.7 

143 

6.4 x 10- 5 

40 

45 

#3 

50 

2.48 

0.27 

100 

+5.2 

-5.2 

3.6 

136.5 

6.8 x 1(T 5 

40 

45 

#4 

50 

2.25 

0.23 

95 

+5.2 

-5.2 

4 

141 

5.5 x 1(T 5 

40 

45 

#5 

60 

2.1 

0.195 

100 

+5.2 

-5.2 

4 

141 

5.5 x 1CT 5 

50 

30 

#6 

80 

2 

0.17 

100 

+5.2 

-5.2 

3.7 

138 

6.4 x 1CT 5 

45 

10 


2 Altitude (km) 



Fig. 6. Normalized electric and magnetic field profiles for the case-study #3 in 
Table 2. Field strength units are for coA 0 = l (Vm 1 ). Amplitudes of horizontal 
electric and magnetic fields are root mean square of 0 and q> components. 



Fig. 7. Comparison between measured ( asterisks and interpolated line ) and retrieved 
electric field strength assumed to be collected by the PWA antenna for three case- 
studies #2, 3 and 4 from Table 2 (dash, dot, and solid lines, respectively ). The retrieved 
profiles are normalized in order to match their peak amplitude with the experi- 
mental value. 


~62 pT. Such surprisingly high values are nearly two orders of 
magnitude larger than the usual values for the Earth’s SR, as 
pointed-out by Simoes et al. (2007). Using the geometric scale 
factor derived in Section 5.1 from the stimulated ELF emissions 
performed by Me Carrick et al. (1990), we saw that a magnetic 
induction of 62 pT would be necessary indeed at the top of the 
cavity to feed the SR mode. Then, the source current-density 


according to the expression )IqI = B would be about 50 pA m ’. 
The latter value is equivalent to a current sheet of 150 A through 
60° in longitude (~2800 km), which is roughly 0.1% of the global 
current induced by the Kronian magnetospheric plasma flow 
(Beghin et al„ 2007). 

5.5. Constrained model of the upper ionosphere conductivity 

Since the atmospheric conductivity profile is seen to control 
Titan’s SR modes, the models built up in Section 5.3 are assumed 
to involve global azimuthal and latitudinal extents, well beyond 
the region of Huygens descent. Therefore, we were forced to re- 
lease some of constraints imposed by the local measurements of 
conductivity by PWA instruments above the GCR layer in favor of 
the constrained parameters derived in Section 5.3 and subse- 
quently introduced in model profiles such as the one plotted in 
Fig. 8 ( solid line). However, although the retrieved field profiles 
plotted in Fig. 7 match satisfactorily the experimental one, we have 
so far no measurement of electron conductivity in the global 
atmosphere up to the upper boundary h 2 to confirm or dispute 
our model profiles. Thus, we were obliged to resort to theoretical 
models. 

Sharp variations of conductivity abovelOO km are thought to be 
linked to haze stratifications distributed in mid latitude regions, 
with scale heights of few kilometers, such as those revealed by 
far-infrared remote sensing from the Cassini CIRS, the Cassini 
Composite Infrared Spectrometer experiment (de Kok et al., 


Altitude (km) 



Fig. 8. Attempts in matching measured and predicted electron conductivity 
profiles, for case-study #5 with parameters from Tables 2 and 3. Asterisks: averaged 
PWA-M1P-RP measurements up to 100 km. Crosses: Model A, with predominating 
GCR ionization, this work. Solid line: model used for retrieval of profiles plotted in 
Fig. 7. Stars: Model B, with predominating solar photoemission threshold of 6.6 eV, 
after Borucki and Whitten (2008). Squares: Model C, same as Model B, but 
photoemission threshold 6 eV. 
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Table 3 

Ionospheric parameters used for the conductivity profiles models A, B and C in Fig. 8. 


Altitude (km) 

Temp. (K) 

Solar attenuation 
(at 250 nm) 

GCR ionization 
rate (cm -3 s -1 ) 

Photoemission threshold (eV) 

Particle size (pm) 

Number density (m 3 ) 


B 

C 

A 

B 

C 

A 

B 

C 

40 

70.66 

5.32 x 1(T 8 

0.45 

n.a. 

n.a 

3.4 

n.a. 

n.a. 

5.0 x 10 6 

n.a. 

n.a. 

50 

71.33 

4.42 x 1(T 7 

0.50 

n.a. 

n.a. 

3.4 

n.a. 

n.a. 

5.0 x 10 6 

n.a. 

n.a. 

60 

76.87 

3.67 x 10- 6 

1.00 

n.a. 

n.a. 

3.4 

n.a. 

n.a. 

5.0 x 10 6 

n.a. 

n.a. 

70 

103.94 

3.04 x 10- 5 

3.02 

n.a. 

n.a. 

3.4 

n.a. 

n.a. 

5.0 x 10 6 

n.a. 

n.a. 

80 

123.61 

2.53 x 10 4 

10.02 

n.a. 

n.a. 

3.4 

n.a. 

n.a. 

5.0 x 10 6 

n.a. 

n.a. 

90 

134.90 

1.79 x 10~ 3 

12.00 

n.a. 

n.a. 

3.4 

n.a. 

n.a. 

4.3 x 10 6 

n.a. 

n.a. 

100 

142.06 

9.62 x 10~ 3 

14.00 

6.6 

n.a. 

3.4 

0.9 

n.a. 

3.7 x 10 6 

1.2 x 10 7 

n.a. 

110 

147.31 

4.06 x 10~ 2 

10.79 

6.6 

n.a. 

n.a. 

0.9 

n.a. 

n.a. 

1.2 x 10 7 

n.a. 

120 

151.00 

1.40 x 10-' 

8.31 

6.6 

6.0 

n.a. 

0.9 

0.9 

n.a. 

1.2 x 10 7 

1.2 x 10 7 

130 

156.34 

4.03 x 10-' 

6.43 

n.a. 

6.0 

n.a. 

n.a. 

0.9 

n.a. 

n.a. 

1.2 x 10 7 

140 

158.55 

9.99 x 10-' 

4.95 

n.a. 

6.0 

n.a. 

n.a. 

0.9 

n.a. 

n.a. 

1.2 x 10 7 

150 

161.92 

1.00 

3.82 

n.a. 

6.0 

n.a. 

n.a. 

0.9 

n.a. 

n.a. 

1.2 x 10 7 

160 

163.77 

1.00 

2.60 

n.a. 

6.0 

n.a. 

n.a. 

0.9 

n.a. 

n.a. 

1.2 x 10 7 


2007). According to modeling of Titan’s atmosphere electric prop- 
erties (Borucki and Whitten, 2008), it appears that the key factors 
controlling the conductivity profile between the GCR layer and 
about 150 km are haze structures of aerosols, whereas little is 
known about size, abundance, photoemission, and electrophyllic 
properties of these particles (Tomasko and West, 2009). Measure- 
ments of the IRIS infrared imager onboard Voyager 1 (Smith et al., 
1981) suggested that highly scattered particles, almost certainly 
condensed organics, with radii between 1 and 5 pm and mole frac- 
tion 10~ 8 -10~ 6 such as nitrile condensates, are presumably perme- 
ating the stratosphere in aerosol clouds at altitudes of 58-90 km 
(Mayo and Samuelson, 2005). However, the most advanced models 
do not account for such thin and irregular structures (Borucki and 
Whitten, 2008). In the latter work, the authors consider particles 
with downward constant mass flux of 4 x 10~ 13 kg nr 2 s~ 3 from 
an altitude of 150 km or above, down to low altitude regions. 
The particle mass density is supposed to be 500 kg m~ 3 , i.e. similar 
to that of tholin aerosols coated with liquid methane (420 kg nr 3 ) 
or ethane (545 kg nr 3 ). 

Stimulated by the encouraging results of Borucki and Whitten 
(2008), predicting electron conductivities at 140 km compatible 
with our constrained values, we tried to identify which parameters 
might be adjusted within superimposed step-by-step structures, in 
order to fit an average profile constrained by the SR model. The re- 
sults are plotted in Fig. 8 with the parameters of the case-study #5 
in Table 2, and three portions (A, B, and C) of profiles with the 
parameters defined in Table 3. The model A, from low altitude up 
to 100 km, appears globally to account quite satisfactorily for the 
PWA-MIP-RP measurements of electron-ion conductivities up to 
the GCR layer (Hamelin et al., 2007; Lopez-Moreno et al., 2008), ex- 
cept the sharp decrease at around 75 km. In model A we have con- 
sidered particles with large aggregate surface radii of 3.4 pm 
(Lavvas et al., 2010), and GCR being the main source of ionization 
(Table 3). Models B and C above 100 km are excerpts from Borucki 
and Whitten (2008) for particles with different photo ionization 
thresholds. Although further refinements could be conceivable, 
the above estimate allows us to check that our constraints applied 
to the upper atmospheric conductivity could be met. We did not 
succeed, however, to simulate the deep decrease of conductivity 
at around 75 km which appears linked to the zonal wind shear 
(Fig. 2), as it falls out of scope of this article, although some related 
aspects are discussed in next Section 6.1. 

6. Discussion 

Among the results of this work, the inferred subsurface struc- 

ture is undeniably the one that is challenging most the models of 

Titan interior. This point is addressed hereunder in Section 6.2, 


while experimental aspects relevant to the specific behavior of 
the PWA experiment are discussed in Appendix B. Concerning 
the atmosphere, we investigate what kind of physical processes 
could possibly account for the apparent correlation between the 
data from independent instruments, as plotted in Fig. 2. Because 
the tilt of Huygens vertical axis with respect to the zenith seems 
to be somewhat correlated with the zonal wind and to have a 
presumed influence on the behavior of PWA boom antenna, we 
discuss the relevant points in Appendix B (Appendix B.3), while 
involved physical processes are analysed here below. 

6.1. Could a correlation between profiles of Huygens tilt, zonal wind 
and conductivity be explained? 

The DWE Team has always been wary of the coincidence (M.K. 
Bird, University of Bonn, private communication) that Huygens 
descended through a surprisingly low-wind region at heights 
100-60 km (Fig. 2 left panel), from a few tens of seconds to several 
minutes after the parachute exchange. In previous studies, it was 
not completely excluded that abnormal parachute dynamics dur- 
ing this phase might possibly be responsible for the apparent drop 
in the zonal wind speed (Bird et al., 2005). On the basis of measure- 
ments of independent instruments, it is proposed here that such an 
artifact is quite unlikely. First, we may rule out any significant link 
between the Huygens vertical axis tilt and the average conductivity 
profile measured by PWA-MIP-RP. Although the mutual impedance 
measurement is influenced to the first-order by the projection of 
the descent velocity vector along the Huygens vertical axis X p 
(Hamelin et al., 2007), second-order effects due to tilt variations 
smaller than 10° would be undetectable. Thus, the correlation be- 
tween the profiles of conductivity and zonal wind is not tilt depen- 
dent. Secondly, we do not see any possibility for the tilt being 
directly caused by the presence of aerosols. 

There is a good correlation, however, between the wind and tilt 
data (Dzierma et al., 2007). If one knows the wind and descent 
speed as a function of altitude, with the probe parameters one 
can predict how much and how long the probe would tilt wherever 
the wind is not constant with altitude. We found that this predic- 
tion is pretty close to the observations under the large parachute, 
at altitudes above 110 km. Under the small parachute, at altitudes 
below 110 km, there is still a small correlation between the 
expected and observed tilt, but there is also another large compo- 
nent of the tilt, which is almost certainly the one-second swing 
when the probe-parachute system oscillates in a scissor-like mode 
(Karkoschka et al., 2007). Although the relatively low value of the 
long time-scale tilt (Q < 5°) was satisfactorily modeled using the 
AGC data from transmission Channel B between Huygens and 
Cassini (Dzierma et al., 2007), the fast oscillations were not 
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expected before the mission and had not been simulated. The latter 
authors suggested that these oscillations could have been 
generated by the probe, and then were picked up by the parachute 
some fraction of a second later, as for instance if one or more para- 
chute risers were kinked by wind gusts or temporarily tangled. 
Wind and tilt are thus somewhat correlated, but not perfectly. 
These observations favor the conclusion that the process responsi- 
ble for the large-scale Doppler wind profile, more especially the 
drop at 75 km, is most likely of atmospheric origin, as it is 
suggested below. 

Concerning the zonal wind profile itself, two reports have been 
subsequently published (Folkner et al., 2006; Flasar et al., 2010) 
after the initial publication of the DWE measurements by Bird 
et al. (2005). Except for an interesting suggestion by Walterscheid 
and Schubert (2006) that atmospheric gravitational tides could be 
responsible for decelerating the superrotating atmosphere at the 
Huygens site, there is still no clear identification of the agent 
responsible for the wind speed drop to nearly zero at 75 km 
(M.K. Bird, University of Bonn, private communication). The mech- 
anism suggested by Walterscheid and Schubert (2006) might ex- 
plain the formation of thin and sharp horizontal haze layers 
assumed to be produced by vertical transport of aerosol particles 
driven by Saturn gravitational tides. According to the standard 
model of superrotating (eastward) zonal wind in equatorial and 
polar regions considered by these authors, the tide-induced verti- 
cal wavelength may become comparable to thin layers of aerosol 
concentration a few km thick, at a critical transition region be- 
tween 90 and 100 km altitude. It is then plausible that similar 
sharp transitions in mid latitude regions could correspond to steep 
gradients of retrieved haze extinction profiles of infrared spectra 
observed by the Cassini CIRS experiment (de Kok et al., 2007). Thus, 
such layering mechanisms, producing compaction and rarefaction 
structures of aerosol concentration in the vicinity of the iono- 
spheric boundaries hi and h 2 , should be considered in further 
development of conductivity profile models, as suggested in Sec- 
tion 5.5. 


6.2. Constrained depth of the buried liquid water ocean and Titan's 
interior models 


The second major constraint derived from our analysis is the 
thickness of the subsurface icy crust (Section 5.2). Assuming the 
lower conductive layer of the cavity to be the crust-ocean inter- 
face, let us summarize the implications of a crust 40-80 km thick. 
Sotin et al. (2009) demonstrated that different models support 
such a thickness range whether internal heat is transferred either 
by conduction or by convection. After describing the consequences 
for each case, we will draw the subsequent inferences on heat flux, 
mantle temperature, and crust composition. If heat is transferred 
by conduction and assuming no lateral variations, the conservation 
of heat can be written 


d_ 

dr 


dT\ 2k dT _ dT „ 

K d7 + T7F = pC ”dt- H - 


2k dT 
r ~dr 


(32) 


where k is the thermal conductivity of the crustal material, r is the 
radial distance, p the density, C p the specific heat, T the temperature 
and H the internal volumetric heating rate. This equation is solved 
with the following assumptions: (i) the thermal conductivity of 
the crust is that of pure ice, (ii) the crust is in thermal equilibrium 
with the internal heat, and (iii) spherical geometry can be neglected 
because the thickness of the crust is small compared to the radius. 
With these conditions, the heat flux through the crust is constant. 
The thermal conductivity of ice strongly depends on temperature 
as k = 0.4685 + 488.1 2/T, which can be approximated at better than 
5% by k = 570/T in the involved temperature range 90-273 K. Note 


that the thermal conductivity varies from k = 6.33 W itT 1 I< 1 at the 
surface (T 0 = 92 1<) to k = 2.1 1 W itT 1 K 1 at the melting tempera- 
ture of pure ice. The 1/T variation leads to the following analytic 
solution of Eq (32): 

in {%) = I in G?) and = 4 ™ 2 ° q = 4na ° in G?) - 

(33) 

where £ = a 0 - r is the depth, z c is the crust thickness (here between 
40 and 80 km), a 0 is the Titan radius, T c the surface temperature, T w 
the temperature of the ocean at the base of the ice crust, q and Qr 
are the heat flux and the total escaping heat power from Titan sur- 
face, respectively. The heat to be transferred through the ice crust is 
the internal heat produced by the decay of radioactive elements 
contained in the rock fraction. This value can be debated, but one 
can take end-members corresponding to the amount in enstatite 
chondrites and carbonaceous chondrites (e.g. Grasset et al., 2000). 
These end-member values provide heat power between 400 GW 
and 650 GW (Fig. 9). If the temperature at the base of the crust is 
the melting temperature of pure ice (T w = 260 I<), then the heat 
which can be released is larger than the maximum value of 
650 GW for values of the crust thickness up to 80 km. Therefore, 
the present values of 40-80 km do not support an ocean tempera- 
ture of 260 K. Note that this value is less than the 273 I< at room 
pressure because the melting temperature of ice decreases with 
increasing pressure. However, the temperature of the ocean can 
be lower, because ammonia, a strong anti-freezing compound, is 
likely to be present in Titan’s ocean (Tobie et al., 2005). With an 
ocean temperature of 200 K, the equilibrium total heat power is 
smaller than 650 GW for values of the crust thickness larger than 
55 km (Fig. 9). The temperature at the crust-ocean interface has 
to be as low as 160 K in order to get a value of 650 GW with a 
40 km thick crust. The constrained range of the crust thickness 
(40-80 km) would be therefore consistent with a crust conduction 
profile over an ammonia-rich ocean with the amount of ammonia 
decreasing (i.e., ocean temperature increasing) with increasing 
thickness. 

The internal heat may also be transferred by convection pro- 
cesses if the buoyant forces overcome the viscous forces. As shown 
by previous studies (e.g. Grasset and Sotin, 1996), convection in Ti- 
tan’s crust can operate when the thickness of the ice layer becomes 
larger than a critical depth (Fig. 10). This critical depth depends on 
the viscosity of ice at the ocean-crust interface: the smaller the vis- 
cosity, the smaller the critical depth. With a typical viscosity of 
10 14 Pa s for ice at its melting point, the critical value of the crust 
thickness is only 15 km. The critical viscosity is equal to 
2 x 10 15 Pas and 2 x 10 16 Pas for a crust thickness 40 km and 

Heat power (GW) 
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Fig. 9. Equilibrium heat power with respect to the crust thickness for two values of 
the temperature at the crust-ocean interface. The thin vertical solid lines are the 
upper and lower bounds of the crust thickness constrained by the SR. The two 
horizontal lines represent end-member values of the internal radioactive power. 
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Fig. 10. Ice viscosity at the crust/ocean interface as a function of the crust thickness. 
The dark grey area corresponds to a domain where heat is transferred by 
conduction (excerpt from Sotin et al., 2009). 

80 km, respectively. Two other elements must be kept in mind. 
First, there are uncertainties on the value of the ice viscosity at a 
given temperature, and the viscosity at the melting temperature 
can vary between 10 13 and 10 15 Pa s. Second, a Titan’s primordial 
large eccentricity implies that there was little dissipation in the 
interior, which yields viscosity larger than 10 15 Pa s at present time 
(Tobie et al., 2005). With such a value of viscosity, a thickness of 
40 km would imply a heat power much larger than 600 GW 
(Fig. 10). On the other hand, a thickness of 80 km is consistent with 
a heat power between 450 and 600 GW and suggests viscosity at 
the base of the crust between 3 x 10 15 and 2 x 10 16 Pa s. The pro- 
posed low bound of 40 km seems therefore difficult to reconcile 
with a convective model, unless additional heat has to be trans- 
ferred such as latent heat released as the ice crust thickens. More- 
over, the upper bound of 80 km is consistent with a convective 
model in equilibrium with radioactive internal heating. Assuming 
that the viscosity law follows an Arrhenius-like law with an 
activation energy of 50 kj/mole and a typical viscosity of 
10 14 Pas at the melting point of pure ice, the viscosity range 
3 x 10 15 -2 x 10 16 Pa s at the crust-ocean interface translates into 
a temperature range 235-220 K, which is in agreement again with 
the presence of ammonia in the ocean. 

The Cassini observations have not yet provided definitive an- 
swers regarding the question about whether convection or conduc- 
tion is the heat transfer mode within the icy crust. For example, the 
presence of 40 Ar in Titan’s atmosphere (Niemann et al., 2005) can- 
not be explained if conduction has operated for several billion 
years. On the other hand, since 40 Ar comes from the decay of 40 I< 
contained in the silicates of the core, there must be some exchange 
between the interior and the atmosphere through convection pro- 
cesses. One may also explain the presence of 40 Ar in the atmo- 
sphere during a global resurfacing event that occurred some 
hundred of million years ago (Tobie et al., 2005). Another fact is 
the absence of evidence for cryo-volcanism at present time. It does 
not prove the absence of convection processes, since other plane- 
tary bodies, like Mars, are likely to undergo convection without 
expressing it by volcanic activity at present time. In conclusion, 
it is not possible today to solve the question of convection versus 
conduction based upon the atmospheric concentration of 40 Ar, nei- 
ther the absence of cry-volcanism on Titan’s surface. 

7. Summary and conclusion 

From an analytic development of mathematical and physical 
properties of the postulated Titan’s SR, in significant progress with 
respect to previous works (Beghin et al., 2007, 2009, 2010), we de- 
rived a model which fits satisfactorily the Huygens PWA-ELF field 


profile measurement. The specificity of Titan’s cavity, with a lay- 
ered inhomogeneous structure and ionospheric current sources in- 
stead of lightning, makes a significant difference with respect to 
the Earth’s traditional SR. Taking advantage of this specificity, in 
spite of severe technical limitations of the PWA experiment, the 
model allows us to constrain the main physical parameters of both 
ionospheric and subsurface boundaries of the cavity. The icy-crust 
properties constrained by this study suggest a quasi-solid sheet 
with finite viscosity overlaying an ammonia-rich ocean, which is 
consistent with convection or conduction processes propagating 
through it. It appears that a crust thickness larger than 80 km 
would be unstable to convection processes unless the viscosity at 
the base of the crust would be larger than 2 x 10 16 Pa s (Fig. 10), 
and the crust-ocean interface temperature less than 260 K. On 
the other hand, a convective 40 km thick crust cannot be in equilib- 
rium with the internal heating, suggesting that the crust should 
continue to grow at the expense of the ocean. Our conclusion is 
that there is, either a ~55 km thick conductive crust with a 
crust-ocean interface temperature of ~200 K and a crust viscosity 
larger than 5 x 10 15 Pas (Fig. 10), or an 80 km thick convective 
crust with an interface temperature of 220-235 K and a maximum 
viscosity of ~2 x 10 16 Pa s (Figs. 9 and 10). Although each of those 
two possibilities of heat transfer would comply with our model of 
the Titan’s SR, it is shown (Section 5.2 and Fig. 4) that a better esti- 
mate or further measurements of both crust permittivity and 
atmospheric conductivity profile, would allow to reduce accord- 
ingly the range of uncertainty. 

After Huygens results and those obtained during the first phase 
of the Cassini mission, we must now wait for long term exploration 
of the Saturn system in order to progress with the question of the 
crust and ocean depth. Unless, in the meantime, the extended Cas- 
sini mission until July 2017 allows advanced measurements of Ti- 
tan's gravity field. It is known indeed that accurate measurements 
of the gravity field might provide the value of the degree 2 Love 
number k 2 , a coefficient of the periodic gravitational potential 
(Sotin et al., 2009). Since this number is very sensitive to the 
presence of a liquid layer at depth, with expected values as large 
as 0.4 compared to values of 1CT 2 when there is no liquid shell 
(Castillo et al., 2000), we may expect from this method, either a 
confirmation that a water ocean is present at depth, or otherwise, 
a challenge to our interpretation of Titan’s SR. 
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Appendix A. Helmholtz equations of Titan’s LSM modes 

While the wave equation is established in its general form (Eq. 
(5)), no special relationship is imposed between both potentials A 
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and *F. Introducing the condition of radial magnetic field compo- 
nent to be zero for TEM modes (Galejs, 1972), and assuming the 
vector potential to have three independent components such as 
A r , A 0 , A v , the radial component of the first equation in Eq. (3) reads 




1 

r sin 0 


d_ 

de 


(A v sin 9) 


d A 
dcp 


(Al) 


Applying the principle of variable separation, it is straightforward 
that the unique solution of Eq. (Al ) implies that both components 
A 0 , and A v are identically null. Thus, the vector potential is radial 
and may be written as a scalar function A=A r (r) Y(0, cp), where 
the spherical harmonic function Y(0, cp) is defined in terms of Legen- 
dre polynomials in Section 5.3. Similarly, the scalar potential reads 
= ‘P r (r) Y(0, cp). Developing each spherical component in Eq. (5), 
yields 


r => 


'V x V x r A r Y(9, cp)' 
r 


,k 2 n 2 d+ r 
co dr 


+ k 2 n 2 A r Y(9, cp) 


0, (A2) 


dr 


- 1 - 


co 


v P r = 0; (p . 


dr 


- 1 


(h 

to 


T, = 0. 


(A3) 


We note that both latitudinal and azimuthal components in Eq. (A3) 
lead to the same relationship between A r and *P r . This relation is 
reminiscent of a gauge condition between both potentials. However 
it is not the conventional Lorentz gauge introduced occasionally as a 
complementary condition (e.g., Sentman, 1990b) that reads as 

k 2 n 2 

V-A-i + r = 0. (A4) 

to v ' 

where the divergence operator replaces the radial derivative in Eq. 
(A3). A gauge condition between both potentials might be arbi- 
trarily chosen as long as no peculiar condition is imposed (Nisbet, 
1955). However, since we imposed radially polarized modes with 
H,-= 0, the “gauge-like” relation given by Eq. (A3) derives directly 
from Eq. (5), and then implicitly from the Maxwell equations as 
shown by Bahar and Fitzwater (1983, Appendix A). Developing 
the first member of Eq. (A2) yields 


VxVx rA r Y(8, cp)' 
r 

_ 1 [ 10/. 0A r Y\ 1 d 2 A r Y 

~ r 2 sine00V Sln ~de~) sin 2 0 d<j> 2 ’ 


(A5) 


which is known as the angular Laplacian of the scalar function A r (r) 
Y{0, cp), and may be developed into spherical harmonics function 
(Groemer, 1996) as 


"V x V x r A r Y(9, cp)' 
r 


1(1+ 1) 

r 2 


A r Y(9, cp), 


(A6) 


where Y(0, cp) is defined by Eq. (30) in Section 5.3. Combining Eqs. 
(5) and (A6) yields a second height-gain differential equation be- 
tween both potentials, in some way symmetrical to Eq. (A.3) 
although independent of it 
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Using both height-gain Eqs. (A3) and (A7) and substituting succes- 
sively Ay then T we obtain two independent second order differen- 
tial Helmholtz equations for the SR Titan’s modes 

i rr+ 2 + (n 2 Vn- 2 ) d ^ = 0. 

k 2 a 2 ) k 2 a 2 n 2 - 1(1+1) dz 


0 2> P r 

dz 2 


+ k l 
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+ k 2 
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k 2 a 2 


A r + (n 2 Vrr 2 ) 


dAy 

Hz 


= 0, 


(A8b) 


in which we made the change of variable z=r - a + z c , where a is the 
mean Titan’s radius and z c the crust thickness (Fig. 1). Moreover, in 
the derivation of these equations we made use of the thin shell 
approximation (Galejs, 1972) which leads to consider k 2 r 2 s» k 2 a 2 
as a constant as far as the mean radius of the cavity is much larger 
than its depth. 


Appendix B. Experimental open questions 

BA. Detection of other eigenmodes 


Here we address the question as to why the PWA instrument 
did not detect any other resonance than the second eigenmode 
at 36 Hz (Section 5.1), although a faint signal above random noise 
of the power-spectrum bin at around 48 Hz (Fig. 7c in Grard et al. 
(2006)), might possibly be a fleeting signature of the mode 1 = 3. A 
plausible cause for the absence of highest order modes could be 
that the Huygens descent latitude (0) would occur for minimum 
values of the spherical harmonic function (Eq. (29)) and its deriva- 
tive, leading to field amplitudes smaller than the background noise. 
Alternatively, a most likely explanation implies that the upper cut- 
off frequency of the ion-acoustic instability (Section 5.1), i.e., the 
plasma frequency of positive ions, would lie at around 40-50 Hz 
in the source region. The ion-plasma frequency obeys the relation 


fi 


e 

2n 



(A9) 


where n± and m, are the positive ions density and mass, respec- 
tively. Assuming that the best conditions for an efficient modulation 
of the ionospheric currents would occur within 200-400 km alti- 
tude range (Beghin et al., 2009), with an electron-ion density be- 
tween 10 6 and 10 7 itT 3 (Borucki et al., 1987, 2006) and a mean 
ion mass between 17 and 170amu, the cut-off frequency would 
be 40-50 Hz, then, ruling out accordingly the source of SR higher 
modes. 

However, the above process cannot explain the absence of the 
fundamental mode (1=1). It is seen from Eq. (25) that the coeffi- 
cient yi is roughly proportional to the ratio z 2 /zi. On Earth, the val- 
ues of the coefficient y, of first two modes are such as 71/72 ~ 1-08. 
As the radial thickness of Titan’s cavity is about twice as large as 
that of Earth, this ratio might even be smaller, say ~1.07. Since 
the coefficient 7 / for the second mode has been prescribed as 
72 = 1.56 (Section 5.2), assuming the ratio 71/72 lying between 
unity and 1.07, the maximum range for 71 is then 1.56 tol.67. 
Substituting these values in Eq. (10), the fundamental mode /1 is 
found between 20.1 and 20.8 Hz. We recall that, due to the loss 
of telemetry channel A (Lebreton et al., 2005), only two spectral 
bins (18 and 24 Hz) among those that might possibly contain /1 
were transmitted to Earth in the highest resolution telemetry 
mode (Fig. 3, left panel). Thus, since /1 lies within the 3 dB band- 
width (21 ±1.5 Hz) of the missing power-spectral bin, its contribu- 
tion could not be observed in the data set, unless overlapping the 
adjacent bins at 18 and 24 Hz significantly above the ambient 
noise. 

It is seen from Eq. (28), that the quality factor for 1 = 1 would be 
roughly the same as for l = 2, yielding a 3.5 Hz half-power band- 
width, that could overlap the adjacent bins with some 3 dB side- 
band damping. Applying to the fundamental mode the same retrie- 
val procedure as for the second harmonic (Section 5.4), keeping the 
same parameters in Table 2 and a constant value for A 0 (Section 
5.1), we find that the fundamental mode would be damped by 
about 6 dB compared to mode l = 2. Adding the noise level of bins 
18-24 Hz (Fig. 3, left panel), the sideband contributions of a 21 Hz 
resonance would be then damped by 6 + 3 dB with respect to its 
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peak value. Thus, such marginal contribution of fundamental Ti- 
tan’s SR could not be seen in the transmitted data. 

B.2. Measurements after landing 

Another challenging issue is why the 36 Hz signal disappeared 
16 s after the Huygens Probe landing. The last significant power 
spectrum (Fig. 3, right panel), was transmitted about 8 s after 
touchdown that occurred at T 0 ~ 8869.80 s Mission Time (Zarnecki 
et al., 2005). The plot in Fig. 3 is an average power spectrum of two 
sample ELF electric field waveforms collected and processed on- 
board about 4 and 6 s respectively after T 0 . Further measurements 
performed 8 s later and transmitted 16.575 s after T 0< and all suc- 
cessive data, exhibit nothing else than instrumental noise which 
reflects a non-nominal performance of the antenna device. This 
fact is corroborated by hard to explain coincidental abnormal data 
of the Mutual Impedance measurements of surface conductivity 
during the first 20 s following touchdown. During the same period, 
three photometers, as a part of D1SR experiment with different 
fields of view, observed an unexplained variability of the light 
intensity for 6 s after T 0 (Chapter 7.1 in Schroder (2007)). 

In short, it seems that Huygens decelerated a fraction of second 
for about 12 cm upon impact while the bottom of the gondola was 
sinking into a soft solid surface (Zarnecki et al., 2005). More pre- 
cisely, the probe bounced back after the first impact, had a free fall 
phase again, then impacted once more and slid for a while, before 
finally coming to rest. The whole sequence lasted slightly less than 
4 s; afterwards Huygens remained motionless (Bettanini et al., 
2008). Thus, the data plotted in Fig. 3 (right panel), were collected 
when the Huygens Probe was at rest with its vertical axis tilted by 
3°, whilst it was only 1-2° during the final part of the descent 
(Karkoschka et al., 2007). The 12 cm penetration may have been 
experienced either by the lower dome of the Probe into the soft 
surface or on ice pebbles that were situated under the dome just 
before touchdown (Bettanini et al., 2008). Therefore, the final 
tilting position might possibly be due to a large icy-cobble on the 
DISR side of Huygens X p axis, like those seen in the HRI and MR1 
images (Karkoschka et al., 2007). It is then plausible that at least 
one of the sensors located at the tip of each boom (Hamelin 
et al., 2007) might have sunk into the ground deeply enough to di- 
vide the initial amplitude of E r by the ground permittivity Re s c as 
seen in Fig. 6. The signal strength would have then decreased by at 
least 6 dB with respect to the records before impact, and hidden 
below the background noise (Fig. 3, right panel). In either case, 
the visible 36 Hz line testifies that the signal was present at Titan’s 
surface at least for 4 s after touchdown. 

B.3. Wind induced vibrations of the PWA boom-antenna 

There was persistent doubt since the earliest investigations that 
the 36 Hz signal could be an artifact somehow related to a mechan- 
ical resonance of the boom antenna (Beghin et al., 2007). First and 
foremost, after static and dynamic cryogenic tests performed in 
2010, on flight-spare models of boom at ESA/ESTEC facilities, it 
was unambiguously concluded that no mechanical resonance 
was observed in the range 36 ± 10 Hz, at all temperatures between 
100 and 300 K. Second, as seen in Fig. 2 (third panel from left), the 
average signal amplitude decays monotonically between 90 and 
55 km, contrary to the zonal wind that drops suddenly from 
50 m s _1 to nearly zero at around 75 km and then grows up to 
35 m s _1 at 60 km (Fig. 2, left panel). Third, after the drogue para- 
chute deployment (110 km altitude), the Huygens Probe descent 
velocity decreases continuously down to the surface, contrary to 
the signal which reaches a peak value at around 90 km and then 
decreases linearly until 60 km. We know that the break of slope 
at 60 km is correlated with the conductivity profile of the GCR 


layer (Beghin et al., 2009). This effect is magnified by a concomi- 
tant break of slope of the Huygens Probe tilt (Fig. 2, right panel), 
that increases accordingly the subsequent contribution of the ra- 
dial electric component (Fig. 7). Thus, although the tilt might be 
somewhat wind dependent (Section 6.1), this does not imply that 
the wind be responsible for a mechanical resonance of the booms. 
Last but not least, since the 36 Hz signal was still present near the 
surface while both wind and descent speeds were weaker than 1 
and 4.6 m s -1 , respectively, as well as a few seconds after landing, 
we can definitely rule out any boom vibration at 36 Hz triggered 
neither by wind-blasts or descent velocity, nor by any reproducible 
process. 
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